function OptFit_Combine(ave_time,upperB,lowerB,x_mean,name,...
    ylabelstr,a,funcIndex,ns_time,outputIndex,plotindex)
% h = figure; %subplot(2,2,plotindex)
% % hold on
% apc3 = 0.775000000000000*0.7;
% apc4 = 0.690000000000000*0.7;
% apc = [0.100000000000000   0.2000000000000   apc3   apc4];

ha = area(ave_time, [lowerB; upperB-lowerB]','FaceColor',[.8 0.8 0.8]);
hold on;
set(ha(1), 'FaceColor', 'none') % this makes the bottom area invisible
set(ha, 'LineStyle', 'none')
pp0 = plot(ave_time,upperB,'--',ave_time,lowerB,'--');
set(pp0,'Color',[0.4 0.4 0.4])


if funcIndex == 2
    % fit to cannonical function
       pp2 = plot(ave_time(1:4:end),x_mean(1:4:end),'ko',ave_time,yd_can(ns_time,a),'r');
       R = corrcoef(x_mean, yd_can(ns_time,a));
       R(1,2)
else
    % fit to linear function
       pp2 = plot(ave_time(1:4:end),x_mean(1:4:end),'ko',ave_time,yd_L(ns_time,a),'r');
       R = corrcoef(x_mean, yd_L(ns_time,a));
       R(1,2)
end


xlim([0 1]);
set(pp0,'MarkerSize', 10,'LineWidth',2);
set(pp2,'MarkerSize', 10,'LineWidth',3);
% set(gca,'position',apc);
xlabel({'Scaled Time'},'Interpreter','LaTex','FontSize',20)
ylabel({ylabelstr},'Interpreter','LaTex','FontSize',20)

set(gca,'FontSize',16)
switch outputIndex
    case 1
    l=legend(pp2,{'$p^H_{COM}$','$p_{COM}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);
    case 2
        l=legend(pp2,{'$m^{H}_{nsCOM}$','$m_{nsCOM}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);
    case 3
        l=legend(pp2,{'$m^{H}_{sCOM} $','$m_{sCOM}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);

    case 4
        ylim([0.6 0.8])
        l=legend(pp2,{'$l^{H}_{sCOM}$','$l_{sCOM}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);
    case 5
        ylim([0.5 0.8])
        l=legend(pp2,{'$l^{H}_{nsCOM} $','$l_{sCOM}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);
    case 6
        l=legend(pp2,{'$\theta^{H}_{nsk}$','$\theta_{nsk}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);
    case 7
        l=legend(pp2,{'$\theta^{H}_{sk}$','$\theta_{sk}^{d}$'},1, ...
   'Location', 'BestOutside', 'Orientation','horizontal',...
   'Interpreter','LaTeX','LineWidth',3,'FontSize',20);
end
% p = get(l,'position');
set(l,'Box','off');
% set(l, 'position', [0 p(2)-.2 1-0.2 p(4)], 'Box','off');
% print(gcf, '-depsc',name)


function ret = yd_can(t,a)
% function 2
ret = exp(-a(4)*t).*(a(1)*cos(a(2)*t)+a(3)*sin(a(2)*t))+a(5);

end

function ret = yd_L(t,a)
    % function 1
        ret = a(1)*t;
end

end
